単純な巡回セールスマン問題

この春，たかしくんは愛媛県庁に入庁した。 これから市町村役場の方々と一緒に仕事をする機会も増えるだろう。 さっそくみんなに顔を覚えてもらうために，愛媛県内にある市区町村役場等の施設にあいさつ回りに行こうと考えた。 それでは，どのようなルートで訪問すればよいだろうか。

1 データの取得

まず，愛媛県の地図データと市町村役場等の施設の位置情報に関するデータをダウンロードすることから始めよう。 国土交通省国土数値情報ダウンロードサイトから，行政区域データと市区町村役場データをダウンロードする。 ファイル名はそれぞれ，N03-20220101_38_GML.zipとP34-14_38_GML.zipである。

ダウンロードしたZIPファイルを解凍し，Rでシェープファイルを読み込む。

library(sf)

ehime_map0 <- sf::read_sf("N03-20220101_38_GML/N03-22_38_220101.shp", options = "ENCODING=CP932")
office_map0 <- sf::read_sf("P34-14_38_GML/P34-14_38.shp", options = "ENCODING=CP932")

このようにしてもよいが，kokudosuuchiを使うとZIPファイルを解凍する必要がなく，便利である。

library(kokudosuuchi)

ehime_map <- kokudosuuchi::readKSJData("N03-20220101_38_GML.zip")
office_map <- kokudosuuchi::readKSJData("P34-14_38_GML.zip") |>
  translateKSJData()

それぞれの方法で読み込まれるデータはまったく同じではないことに注意が必要である。 後者はリストのリストとして読み込まれる。

is.list(ehime_map0[[1]])

[1] FALSE

is.list(ehime_map[[1]])

[1] TRUE

市町村名をfactorにしておくと，図を描いたときの凡例の順序がおかしくならずにすむ。

ehime_map[[1]]$N03_004 <- factor(ehime_map[[1]]$N03_004, levels = unique(ehime_map[[1]]$N03_004))

1.1 データの修正

先ほどダウンロードした愛媛の令和4年の行政区域データには，市区町村名（N03_004の列）がNAとなっているレコードが存在する。

table(is.na(ehime_map[[1]]$N03_004))

FALSE  TRUE 
 4906     1 

この部分を残したままだと，以下の図の凡例にNAという名前の市町村が現れる。 これをどう処理するかは考える必要があるが，たかしくんは強引にもレコードごと削除した。

ehime_map[[1]] <- ehime_map[[1]][-which(is.na(ehime_map[[1]]$N03_004)), ]

2 地図のプロット

Rにデータを読み込めたので，さっそく地図をプロットしてみる。

# plot(ehime_map[[1]])# 処理に非常に時間がかかるため，避けた方がよい
# plot(sf::st_geometry(ehime_map[[1]]))
library(rmapshaper)

ehime_map_shrink = rmapshaper::ms_simplify(ehime_map[[1]], keep = 0.001, keep_shapes = TRUE)

par(mar=c(0, 0, 0, 0)) 
plot(sf::st_geometry(ehime_map_shrink))
plot(sf::st_geometry(office_map[[1]]), pch = 16, col = "orange", add = TRUE)

ぜんぜん悪くない。 ただし，たかしくんはggplot2を用いた図の方が好みである。

library(ggplot2)
# カラーパレット
library(viridis)

ggplot(ehime_map[[1]]) +
  geom_sf(aes(fill = N03_004)) +
  geom_sf(data = office_map[[1]], colour = "orange") +
  scale_fill_viridis_d(option = "mako", direction = 1) + 
  guides(fill = guide_legend(title = "市町村名")) +
  theme_void() +
  theme(text = element_text(family = "HiraKakuProN-W3"))

たかしくんはここで重大な問題に気がついた。 この地図は動かせないし，拡大もできない。 Google マップに慣れ親しんでいるたかしくんにとって，地図は動かせた方が便利である。 そこで，Leafletを使うことにする。

library(leaflet)

micanfIcon <- makeIcon(
  iconUrl = "https://www.pref.ehime.jp/h12200/mican-kanzume/images/img_dance02.png",
# iconUrl = "https://www.pref.ehime.jp/h12200/mican-kanzume/images/img_dance03.png",
  iconWidth = 40, iconHeight = 40,
  iconAnchorX = 20, iconAnchorY = 20,
)
leaflet::leaflet(office_map[[1]]) |>
  addTiles() |>
  addMarkers(icon = micanfIcon)

みきゃんだ。

library(mapview)

mapview::mapview(ehime_map[[1]], zcol = "N03_004", layer.name = "市町村名") +
  mapview(office_map[[1]], col.regions = "orange", legend = FALSE)

訪問場所が多いので，何とか減らせないかとたかしくんは考え始めた。 松山市内だけにするか，あるいは，市町村役場だけにするか。 前者だとすぐに回れるが，後者の方がいろいろ行けて楽しそうだ。 たかしくんは，両方のルートに必要な移動時間と移動距離を計算してから，どちらにするか考えることにした。

3 ルートの探索

3.1 松山市内だけを回るルート

library(osrm)

matsuyama <- office_map[[1]][grep("松山市", office_map[[1]]$所在地), ]
matsuyama_tsp <- osrm::osrmTrip(loc = matsuyama, returnclass = "sf", osrm.profile = "foot")

mapview::mapview(matsuyama, legend = FALSE) +
  mapview(matsuyama_tsp[[1]]$trip, color = "black", homebutton = FALSE, legend = FALSE)

ルート全体の移動時間は2236.4分，移動距離は201.2kmであった。

# Graphviz（ここでは図示しない。）
library(DiagrammeR)

nodes_df <- DiagrammeR::create_node_df(
  n = nrow(matsuyama),
  label = matsuyama$名称,
  shape = "plaintext",
  color = "black",
  fillcolor = "white")
edges_df <- DiagrammeR::create_edge_df(
  from = matsuyama_tsp[[1]]$trip$start,
  to = matsuyama_tsp[[1]]$trip$end,
  color = "black")
edges_df$label <- paste(
  paste("duration", round(matsuyama_tsp[[1]]$trip$duration, 1), sep = ": "),
  paste("distance", round(matsuyama_tsp[[1]]$trip$distance, 1), sep = ": "),
  sep = "\n")
g <- DiagrammeR::create_graph(
  nodes_df = nodes_df,
  edges_df = edges_df,
  attr_theme = "lr")
DiagrammeR::render_graph(g)

library(visNetwork)

nodes <- data.frame(
  id = rownames(matsuyama),
  label = matsuyama$名称,
  shape = "text"
)
edges <- data.frame(
  from = matsuyama_tsp[[1]]$trip$start,
  to = matsuyama_tsp[[1]]$trip$end,
  label = paste(
    paste("時間", sprintf("%.1f", round(matsuyama_tsp[[1]]$trip$duration, 1)), sep = ": "),
    paste0("距離", ": ", sprintf("%.1f", round(matsuyama_tsp[[1]]$trip$distance, 1)), " "),
    sep = "\n"
  ),
  color = "black", arrows = "to",
  shadow = TRUE
)
visNetwork::visNetwork(nodes, edges) |>
  visOptions(highlightNearest = TRUE, nodesIdSelection = TRUE)

（この図は拡大できます。）

3.2 市町村役場だけを回るルート

cityhall <- office_map[[1]][grep("役所$|役場$", office_map[[1]]$名称), ]
cityhall_tsp <- osrm::osrmTrip(loc = cityhall, returnclass = "sf", osrm.profile = "car")

mapview::mapview(
  list(cityhall, cityhall_tsp[[1]]$trip),
  color = list("black", "black"),
  col.regions = list("red", NA),
  legend = list(FALSE, FALSE),
  homebutton = list(TRUE, FALSE)
)

ルート全体の移動時間は539.5分，移動距離は649kmであった。

nodes <- data.frame(
  id = rownames(cityhall),
  label = cityhall$名称,
  shape = "text"
)
edges <- data.frame(
  from = cityhall_tsp[[1]]$trip$start,
  to = cityhall_tsp[[1]]$trip$end,
  label = paste(
    paste("時間", sprintf("%.1f", round(cityhall_tsp[[1]]$trip$duration, 1)), sep = ": "),
    paste("距離", sprintf("%.1f", round(cityhall_tsp[[1]]$trip$distance, 1)), sep = ": "),
    sep = "\n"
  ),
  color = "black", arrows = "to",
  shadow = TRUE
)
visNetwork::visNetwork(nodes, edges) |>
  visOptions(highlightNearest = TRUE, nodesIdSelection = TRUE)

3.3 ルートの比較

これら2つのルートを比較しながら，たかしくんは考え込んでしまった。 移動距離が短い方が移動時間が長い？ どうしてこうなったか考えてみよう。

	移動時間	移動距離	訪問場所
松山市内だけを回るルート	2236.4 分 （= 37.3 時間）	201.2 km	30 か所
市町村役場だけを回るルート	539.5 分 （= 9.0時間）	649.0 km	20 か所